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Abstract 

Multi-loop integrals can be evaluated numerically using Mellin-Barnes representa- 
tions. Here this technique is applied to the calculation of electroweak two- loop correc- 
tion with closed fermion loops for two observables: the effective weak mixing angle for 
bottom quarks, sin 2 9^, and the branching ratio of the Z boson into bottom quarks, 
Rb. Good agreement with a previous result for sin 2 9^ is found. The result for R), is 
new and it turns out that the electroweak two-loop corrections are relatively large - 
of the same order as the current experimental uncertainty. A simple parametrization 
formula is provided which approximates the full result within integration errors. 



1 Introduction 



High-precision data from LEP, SLC, the Tevatron, and the LHC allow us to perform very 
accurate tests of the Standard Model (SM), predict the mass of the Higgs boson, and possibly 
identify hints for new physics. For this purpose, the experimental results are compared to 
theoretical computations including higher-order radiative corrections. For some of the most 
precisely determined observables it is necessary to include complete two-loop corrections 
in the calculation. Such calculations have been carried out for the mass of the W boson, 
M\v [IH3], and the effective weak mixing angle sin 2 9 e g parametrizing the ratio of vector 
and axial- vector couplings of the Z boson to leptons [IH6], light quarks and bottom 
quarks [7]. Furthermore, partial electroweak three-loop contributions [8], as well as three- 
loop PUIEEI] and leading four-loop [T2] QCD corrections to these observables have been 
determined. The two-loop and leading three-loop results for M w and sin 2 Q eS have been 
implemented in commonly used SM fit programs such as ZFitter [13J and GFitter |14j . 

A major difficulty in these calculations are two-loop integrals with multiple mass and 
momentum scales. In general these integrals cannot be solved analytically in closed form, 
so that one has to use numerical methods instead. Numerical techniques for the evaluation 
of loop integrals face two main challenges: extraction of ultraviolet (UV) and infrared (IR) 
singularities, as well as stability and efficiency of the numerical integrations. A powerful 
method is based on Mellin-Barnes (MB) representations. It has been demonstrated that MB 
representations can be used to isolate the singularities of an arbitrary multi-loop integral in 
a systematic way [To^fTB], such that this procedure can be automated in a computer program 
[TfijfTT]. Furthermore, it was shown that the integration time and convergence behavior of 
the MB integrals can be improved substantially by suitable variable transformations and by 
analytically integrating over some variables of the multi-dimensional MB integrals [T8] . 

This article reports on the concrete calculation of two electroweak precision observables 
using the technique described in Ref. [IE]. First, the method is validated by reproducing 
the result for the fermionic electroweak two-loop corrections to the effective weak mixing 
angle of bottom quarks, sin 2 6^, published earlier in Ref. [TJ. Here the term "fermionic" 
refers to diagrams with at least one closed fermion loop. Secondly, a new result is presented 
for the complete fermionic electroweak two-loop (next-to-next-leading order) corrections to 
Rb, the branching ratio of the Z boson into bottom quarks and all quarks. Until now, only 
an approximate results for the electroweak two-loop corrections to the Z partial widths are 
available, using an expansion for large values of the top-quark mass. For decays of the Z 
boson into quarks of the first two generations, this expansion has been driven to the order 
0(a 2 m 2 ) [19], while for the Z bb width only the leading 0{a 2 mf) coefficient is known [20J. 

As elaborated in the next section, the radiative corrections for both sin 2 9^ and are 
obtained from the same loop diagrams of the Zff vertex. For the calculation presented here, 
the MB method of Ref. [TS] has been used for the most difficult two-loop diagrams involving 
triangle sub-loops. The remaining diagrams have been computed by reducing the relevant 
integrals to a set of master integrals [21] , which are then evaluated numerically [5l|22] . These 
steps are described in more detail in section |3j Finally, numerical results for both observables 
are presented in section HJ and the main findings are summarized in section 
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Figure 1: Fermionic electroweak two-loop Zff vertex diagrams contributing to the vertex 
form factors vf and af . 

2 sin 2 an d R b a t next-to-next-to-leading order 

The effective weak mixing angle sin 2 9^ is related to the ratio of vector and axial- vector 
couplings, Vf and a { of the Zff vertex, i. e. 



4 V a f 



Expanding this definition to next-to-next-to-leading order (NNLO) one obtains 

M 2 \ ( n (1) v {0) - v W n {0) 
sm 2 6{l = ( 1 - Ue 1 + Gf "J ^ Vf af 



k 2 =Ml 



^ af vf qf - v f (aff - {aff vf + af vf af 



(afnaf-vf) 



(2) 



fc 2 =M| . 



where v f 1 ' 2 ^ denote the tree-level, one-loop, and two-loop contributions to the vector form 
factor, respectively, and af 1 ' 2 ^ are defined similarly for the axial- vector form factor. Figured] 

(2) 

shows the types of two-loop diagrams with closed fermions loops that contribute to v) and 

(2) 2 2 

a j . M w and M z are the real parts of the gauge- invariant propagator poles of the W and 
Z boson, respectively. In eq. fl2]) the infrared (IR) divergent QED and QCD contributions 
cancel exactly, see for example Ref. [5], so that it is sufficient to only include IR- finite weak 
loop corrections. 
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The branching ratio Rb is defined as the ratio of the partial decay widths of the Z-boson 
decay into bottom quarks and into all quarks: 

b r had v d + v u + v s + r c + T b 1 + 2(v d + rj/ry 1 ; 

where Tf stands for the partial decay width into the // final state. In the last step in fl3]), the 
relationships T u T c and ~ T s have been used, which hold to very good approximation. 
Up to next-to- next-to-leading order (q = u,d), 



(GH 2 J (G^) 2 

+ -i-^e{ (Gf) 2 [2Gf) + (4 1 )) 2 + (^)) 2 ] 
1 6 j -GfGf)[2G« + (a«) 2 + (t,«) 2 ] 



(1)\2 



+ TT^b [^^f - G W - + Gf(^) 2 ] (4) 

(Gh T 
2 

(g] 

+ G(°)G«^ + G<° W ( W + WX^i - <v) 
-G^G^], 



with 



G (n) = a (0) Q (n) + 17 (0, W ( B)) ^(n) = {a ff n ^ + (,f ) 2 ^. (5) 

Here 7^ n y and incorporate the n-loop QED and QCD corrections to the vector and 
axial- vector form factors, which have been calculated already several years ago [231124] . see 
also Ref. [25J. The relevant parts for this calculation are given by 

jjW _ -r^ 1 ) _ jyi 1 ) — _a_ _|_ * (a) 

/C «,V - /C «,A - 37T + TT ' \' ) 

<! = <v - 6^, (8) 

<v = -^ + G 2 (t) 2 + G 3 (f) 3 , (9) 

^-^ + '.(#) + («', (10) 

tt«s 

12tt 2 ' 

<i=<v-^@(t) z -/3(^)(tr, (12) 

^ 2 v = + 12|§t - 6 S| + K«s, m 2 a 2 ), (13) 



'T->( 2 ) _ "T?( 2 ) acts (M\ 
'VV — /C d,V ~~ 12^' ^ iJ -^ 



<i = <i - 22^f + 6^ + OK« Sj mfc 2 ). (14) 
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Here three-loop QCD corrections are also included in the VS 2 ^ terms, and is defined 
in the MS-scheme. In the power-counting the small bottom quark mass is treated to be 
parametrically of the same order as the fine-structure constant, jp ~ a, so that only the 

leading tree-level mb-dependence is included in ([5} and the leading one-loop mb-dependence 
in (ITS]) . flHJ). Note that the dependence on C 2 and C 3 drops out in the final result for R b . 
The functions I 2 and J 3 read 

37 7 

h{x) = -r^ + log x+—x + 0.0132 x 2 , (15) 
12 ol 

5075 23 . . 67, 23 , , 

J 3 (a;) = 1 C2 + C3H logaH log 2 x. (16) 

; 216 6 s s 18 B 12 B v ; 

(2) (2) 

The vertex form factors v j and include all electroweak two-loop diagrams, as shown 
in Fig. [TJ except those that involve IR-divergent photon exchange contributions between 
the outgoing fermion lines (represented by Fig. [1] (a) when both wiggly lines are photons). 
For the effective weak mixing angle sin 2 9 e s these QED diagrams cancel in the final result, 
as pointed out in Ref. [5]. For the branching fraction Rt,, on the other hand, the QED 
contributions are included in the radiator functions H q y and lZ q ^. In consequence, for 
both observables, the IR-divergent QED diagrams are excluded from the calculation of the 
electroweak two-loop corrections discussed here. 

For the renormalization the on-shell scheme is employed. In particular, the renormalized 
squared masses are defined as the real part of the propagator poles (which is different than 
the pole of the real part of the propagator). Furthermore, the external fields are renormalized 
to unity at the position of the poles. Details and explicit expressions for the renormalization 
constants can be found in Ref. [2]. 

For the computation of the electroweak two-loop corrections, the masses and Yukawa 
couplings of all fermions except the top quark can be safely neglected. Moreover, the CKM 
matrix is approximated by the unit matrix. 



3 Outline of the computation of two-loop diagrams 



Let us start by giving a brief overview to the use of MB representations for the numerical 
evaluation of multi-loop integrals. For more details, see Refs. [T5lfT6l fT8] . 

After introduction of Feynman parameters, a general one-loop integral with iV propaga- 
tors can be written in the form 



l) N Y{N-D/2) I dx x ---dx 



'N 



5(1 — xi 



x N , 



Yaj=i K ij X i X j + J2iLl L i X i + M -ie 

where D is the number of space-time dimensions, and Kij, Li, and M depend on the masses 



^A' 



N-D/2 



, (17) 
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and external momenta of the propagators. Using the MB representation, 

d Zl --- [ dz rn A?---A z -A z -^--: 



(A + ... + A m ) z (2m) m J Ci 



T(-zi) ■ ■ ■ Y(-z m )Y(Z + zi + ...+z n 



Y(Z) 

the integral / can be transformed into a form that depends on the Feynman parameters only 
in terms of exponentials x Zi . The integration contours Ci for Zi are straight lines parallel to 
the imaginary axis chosen such that all arguments of the gamma functions have positive real 
parts. MB representations for multi-loop integrals can be obtained recursively. After inte- 
gration over the Feynman parameters, the remaining MB integrals only depend on gamma 
functions of the integration variables and on the external parameters. For example, the 
"sunset" two-loop diagram is given by the MB integral 

dz x dz 2 dz 3 {ml ) ~ £ - Zl - Z2 {mlY 2 {m 2 3 ) ~ 23 ( -p 2 ) z * 

x r(-z2)r(-z 3 )r(i + z x + Z2 )y(z 3 - Zl ) (19) 

r(l - e - z 2 )Y{e + Zl + z 2 )Y{e -l- Zl + z 3 ) 




(2ttz) 



X 



T(2-e + z 3 ) 

with e = (4 — D)/2. All UV and IR singularities are now in the poles of the gamma functions 
and can be extracted in a systematic way, see Refs. [15l[16] . 

The remaining integrals in Z{ over the complex contours Ci can, in principle, be carried out 
numerically. However, in practice, the integrand is often highly oscillatory and the integral, 
while formally existent, may converge too slowly. 

In Ref. [18], two effective methods for improving the convergence behavior of the numer- 
ical Zi integrals have been discussed. One modification consists in rotating the ^-integration 
contours in the complex plane, which can produce an additional exponentially damped term 
in the integrand, so that the integral vanishes faster for \zi\ — > oo. By rotating the contours 
for all variables Zi in parallel it is ensured that no poles of the gamma functions cross the con- 
tour. Secondly, some integrations of the multi-dimensional MB integral can be performed 
analytically with the help of the convolution theorem for Mellin transforms. For typical 
two-loop vertex diagrams needed for the calculation of the NNLO corrections to the Zbb 
interaction, one can treat about half of the integrations in this way, and only the remaining 
half will be carried out numerically. With these improvements, it was shown in Ref. [18J that 
the most complicated scalar integrals needed for the NNLO Zbb corrections can be evaluated 
in a few hours on a single-core computer. 

In this work, the MB method has been used for the most complicated integrals stemming 
from the diagrams with triangle sub-loops, see Fig. [1] (c). The remaining two- loop integrals 
have been calculated by first reducing them to a set of master integrals using integration-by- 
parts and Lorentz identity relations [21] . Analytical formulas are known for the one-loop and 
two-loop vacuum master integrals, while the two-loop selfenergy and vertex master integrals 
have been evaluated numerically using dispersion relations as described in Ref 
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Parameter Experimental value 



M z 


(91.1876 ±0.0021) GeV 


Tw 


(2.4952 ± 0.0023) GeV 


M w 


(80.399 ± 0.023) GeV 


T\v 


(2.085 ± 0.042) GeV 


m t 


(173.2 ±0.9) GeV 


Aa(M z ) 


0.05900 ± 0.00033 


m h 


4.20 GeV 


«s(M z ) 


0.1184 ±0.0007 




1.16637 x 10" 5 GeV~ 2 



Table 1: Input parameters and their experimental values, from Refs. {2 ty l2T I . The elec- 
troweak two-loop corrections of the Zbb vertex depend only the first group of parameters. 



4 Results 

The computational method described in the previous section has been applied to the cal- 
culation of two- loop corrections to sin 2 9^ and R^. Both quantities depend on the input 
parameters listed in the upper part of Tab. [TJ Note that the experimentally quoted val- 
ues for the W- and Z-boson masses correspond to a Breit-Wigner parametrization with a 
energy-dependent width, and they have to be translation to the pole-mass scheme used in 
the loop calculation [28]. In effect, this translation results in a downward shift of M\y and 
Mz by r 2 v /(2Mw) and r|/(2Mz), respectively, where Tw,z and width of the gauge bosons. 
Aa(Mz) is the contribution from light fermions (all fermions except the top quark) to the 
running of the electromagnetic coupling between the scales Q = and M%. 

Let us first discuss of the effective weak mixing angle of bottom quarks, sin 2 6^. The 
fermionic two- loop corrections to this quantity have already been calculated earlier [7], so 
that this result can be used as an additional check of the MB method. The comparison is 
shown in Tab. [2] in terms of the two-loop correction factor in the second line of eq. (J2J): 



(2) (o) (o) (2), (o)n 2 /_(i)\2„(o) , mm JO) 

A/ > 2 ,ferm) _ % V b a b - V b (d b ) ~ (d b ) V b + d b V b d b 



(20) 



The two calculations agree to within about 1% of the two-loop contribution, which is the 
same order as the uncertainty from the numerical MB integration. 

Next, we turn to the presentation of results for the branching ratio Rb- Table |3] lists 
the effects of various radiative correction terms, starting with the one-loop contribution, 
corresponding to the first line of eq. (T4|), together with C(a 2 ) final-state corrections, which are 
numerically of similar order. The third column of Tab. |3] shows the new result for the purely 
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Mw m t A«jf ,ferm) x 10 3 



[GeV] [GeV] Ref. [7] this work 



80.399 


171.2 


-2.278 


-2.303 


80.399 


172.2 


-2.331 


-2.357 


80.399 


173.2 


-2.386 


-2.410 


80.399 


174.2 


-2.441 


-2.464 


80.399 


175.2 


-2.496 


-2.518 


80.422 


173.2 


-2.402 


-2.427 


80.445 


173.2 


-2.418 


-2.444 



Table 2: Results for the two-loop correction factor Ak^ ' e to sin 2 9^, for different values 
of Mw and m t . The other input values have been set to Mz = 91.1876 GeV, Mr = 100 GeV, 
Aa = 0. The values obtained in this work are compared to the results of Ref \7§. 

electroweak contributions from two-loop diagrams with closed fermion loops. The fourth 
includes the electroweak two-loop corrections together with final-state corrections of order 
0(a 3 , aa s , m 2 a s , m 4 ). Finally, higher-order two-loop QCD corrections of order 0{aa s ) [29] 
to internal gauge-boson selfenergies, as well as three-loop QCD corrections of order (9(aa 2 ) 
to the p-parameter [TU] are given in the last column of the tabled. Additional higher-order 
corrections to the p-parameter [3[9l[lTl[T2] are very small and have not been included in the 
numerical analysis. 

As evident from the table, the electroweak two-loop corrections are relatively large com- 
pared to the one-loop results, and they are comparable to the experimental uncertainty for 
Rb of 6.6 x 10~ 4 [26] . The two-loop contribution depends only mildly on the Higgs-boson 
mass. Its overall effect is a reduction of the SM prediction for R b . It is interesting to note 
that the electroweak two-loop corrections, the final-state QED and QCD corrections, and the 
higher-order QCD contributions to internal gauge-boson lines all produce a negative shift 
of the predicted values for R b . The uncertainty in the prediction for R b due to the numer- 
ical MB integration error is about 4 x 10 -6 , i. e. less than 1% of the electroweak two-loop 
contribution. 

The numerical results in Tab. [3] have been obtained with a fixed value of the H^-boson 
mass. In global SM fits, however, Mw is calculated from the Fermi constant G M . In Refs. p~|- 
[3] this calculation has been carried out including full two-loop and partial higher-order 
corrections. Using the same order of perturbation theory for the calculation of Mw and R b , 
the numerical results for R b in this scheme are given in Tab. H] and Fig. [21 As evident from the 
figure, the higher-order corrections beyond the one-loop level are sizeable. The most precise 
result including electroweak two-loop and QCD three-loop corrections has a mild tension of 
about two standard deviations with the experimental value R b = 0.21629 ± 0.00066 [26] . 

The evaluation of the complete fermionic two-loop result with the MB integrals takes 

*Non-fact ??? 
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Mil 
[GeV] 


0(a) + FSR^ ^2 

[io- 3 ] 


0(a 2 ) 

^ Inform/ 

[io- 4 ] 


0{a} ) + FSR 3 2 4 
[IO" 4 ] 


OiaOc am 2 ) 
[IO" 4 ] ' 


1f)f) 


—3 566 

•J . WWW 


—6 583 

w . WOW 


—8.214 


—0 404 

W . Tfc Wt: 


200 


-3.585 


-6.587 


-8.216 


-0.404 


400 


-3.609 


-6.595 


-8.220 


-0.404 


600 


-3.624 


-6.594 


-8.216 


-0.404 


1000 


-3.645 


-6.582 


-8.201 


-0.404 



Table 3: Results for electroweak one- and two-loop corrections to R b , as defined in eqs. ( BE]) , 
for different values of Mh- Tie other input values are taken from Tab. (TJ with a fixed value 
for M\y. Also shown are the effects of two- and three-loop QCD corrections to the final state 
(fourth column) and to gauge-boson selfenergies (fifth column). Here "FSR" stands for the 
final-state radiative QCD and QED corrections described by the radiator functions IZ^ . 



M H 
[GeV] 


tree- level + 0(a) 
-\- FSR Q Q , sjQ ,2 


0(^ferm) 'PSR ^ > exas,m?a s ,'rni 

+ 0(aa s , aa 2 ) 

[io- 4 ] 


total 


100 


0.21569 


-7.616 


0.21493 


200 


0.21570 


-7.583 


0.21494 


400 


0.21572 


-7.548 


0.21496 


600 


0.21573 


-7.531 


0.21497 


1000 


0.21574 


-7.517 


0.21499 



Table 4: Results for R b , as in Table but now with M w calculated from using the SM 
prediction. The other input values are taken from Tab. [IJ 



several CPU-hours for one set of parameters and thus is not suitable for direct incorporation 
in SM fit programs. However, the prediction for R b can be approximated to a good precision 
by a simple parametrization formula: 



with 



R b = R° b + c x L H + c 2 L 2 H + c z L A H + c 4 (A 2 H - 1) + c 5 A a 
+ c 6 At + c 7 A t L H + c 8 A as + c 9 A 2 + c w A z , 



, M H M H ( m t 

m = In — — — — . A H = A t 



(21) 



100 GeV 100 GeV V 173.2 GeV 

A a = ^_-1, A a =^l-1, A Z = ^^ 1. (22) 

0.05900 ' 3 0.1184 ' 91.1876 GeV V ' 

The numerical coefficients are determined by a fit to the full numerical result, which includes 
all radiative corrections mentioned above: the complete O(o) and fermionic 0(a 2 ) contri- 
butions to the Zff vertex form factors, as well as virtual 0(aa s ) and 0(aa 2 ) corrections 
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0.2158 
0.2156 
0.2154 
0.2152 
0.2150 
0.2148 



200 



1-loop [incl. 0(a, a s , a 2 s ) FSR] 



+ 0(aj erm , aa s , aa s ) [incl. higher-order FSR] 



400 600 
M H [GeV] 



800 



1000 



Figure 2: One-loop and two-loop (with QCD three-loop contributions) result for the branch- 
ing fraction Rb as a function of Mh, and using the SM prediction for M\y to the same order 
of perturbation theory. Input values for the other parameters are taken from Tab. QJ 



and final-state radiation of order 0(o>™), (n = 1,2,3) and 0(aa s ). For the PF-boson mass 
the currently most precise result of Ref. [3] is used. With these ingredients, the fit result for 
the coefficients is 

R° b = 0.2149246, d = 2.23 x 10" 5 , c 2 = 2.6 x 1(T 6 , c 3 = -6.8 x 1(T 7 , 

c 4 = 9.19 x 10" 8 , c 5 = 6.58 x 10" 4 , c 6 = -3.363 x 1(T 3 , c 7 = 6.74 x 10" 5 , (23) 
c 8 = -1.688 x 1(T 3 , c 9 = -9.26 x 1(T 4 , c w = 5.93 x 1(T 2 . 

With this parametrization the full result is approximated to better than 10~ 6 for 10 GeV < 
Mh < 1 TeV and the other input parameters in their la experimental error ranges. 



5 Summary 

Mellin-Barnes representations provide a generic framework for evaluating multi-loop integrals 
numerically. They can used to systematically extract all physical singularities of a given loop 
diagram. Furthermore, by using appropriate variable mappings and contour deformations 
one can achieve good numerical convergence of the Mellin-Barnes integrals for a large class 
of loop diagrams, see Ref. [18] . 

This paper reports on the application of this method to the calculation of fermionic elec- 
troweak two-loop corrections to the effective weak mixing angle for bottom quarks, sin 2 9^, 
and the branching ratio of the Z boson into bottom quarks, Rb- These two observables have 
been measured to high accuracy by experiments at LEP and SLC, and they play an impor- 
tant role in precision tests of the Standard Model, as well as new physics models. Therefore 
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one needs theoretical predictions for these quantities with an error that is comparable or 
less than the experimental uncertainty, which requires the inclusion of electroweak two-loop 
corrections. 

In the work presented here, the numerical Mellin-Barnes method has been used for the 
most difficult two-loop vertex diagrams. For the one-loop and remaining two-loop diagrams, 
a technique based on tensor reduction to a set of well-known master integrals has been 
employed, which is faster but less general than the Mellin-Barnes method. 

The results for sin 2 9^ have been compared to a previous calculation of the fermionic 
two-loop corrections with a different method Ref. [7]. Good agreement within numerical 
integration errors has been found, which serves as a validation of the numerical MB method. 

The result for Rb presented here is new. The fermionic electroweak two-loop corrections 
are found to be relatively large, of the same order as the experimental error of 6.6 x 10~ 4 . 
Thus it is important to include these corrections in global electroweak precision fits. For easy 
use by other groups, a simple parametrization formula has been provided, which accurately 
approximates the full result within integration errors. 
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